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ABSTRACT 

This paper presents a general approach to linear stochastic processes driven by various 
random noises. Mathematically, such processes are described by linear stochastic differen- 
tial equations of arbitrary order (the simplest non-trivial example is x = R(t), where R(t) 
is not a Gaussian white noise). The stochastic process is discretized into n time-steps, 
all possible realizations are summed up and the continuum limit is taken. This proce- 
dure often yields closed form formulas for the joint probability distributions. Completely 
worked out examples include all Gaussian random forces and a large class of Markovian 
(non- Gaussian) forces. This approach is also useful for deriving Fokker-Planck equations 
for the probability distribution functions. This is worked out for Gaussian noises and for 
the Markovian dichotomous noise. 

P. A. C. S. Numbers: 05.40,+j., 02.50.-r. 
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I. Introduction 

Stochastic differential equations were first studied in the context of Brownian motion 
by Langevin. For a particle of unit mass, he derived the equation for the displacement 

x + ~fx = R(t). (1.1) 

In this equation, the drag force, — jx, is a deterministic term which represents the "envi- 
ronment" in which the particle moves. In the case of Brownian motion, this is the fluid 
surrounding the particle. The force driving the particle through this environment is the 
random term R(t) which is assumed to be a Gaussian white noise [1]. 

In recent years, the field of stochastic equations has expanded considerably along two 
main lines, each focusing on one of the two basic aspects of Langevin equations: The first 
concentrates on more complex environments, usually introduced through potentials U(x); 
the second concentrates on random driving forces more complex than the Gaussian white 
noise. 

The first approach deals with equations like, e. g., 

x + ^M=fl(f), (1.2) 

where U(x) is an external potential [2]. Since Eq. (1.2) is typically non-linear and very 
hard to solve, the random noise is usually kept Gaussian and white. (As an aside, note 
that when one tries to combine non-trivial potential terms with colored noises [3,4,5], the 
difficulties are so great that the order of the equation has to be reduced to the form [4] 

x + n(x) = R(t) ). 

The present work, however, follows the second line of generalizations of Eq. (1.1), 
which leaves the "environment" relatively simple but considers more complex random 
forces without reducing the order of the equation. It has been known for some time, 
for example, that a more realistic description of Brownian motion must include a finite 
correlation time for the random force [6], a feature absent in the Gaussian white noise. If 
the random force arises from processes internal to the system (as in Brownian motion), 
the system may reach a state of detailed balance at equilibrium. In this case, if the force 
is not Gaussian white, the friction term requires a modification which usually involves a 
retardation effect [7]. Eq. (1.1) is then generalized to be : 

x+ I 7 (£ - t')x(t') dt' = R(t) , (1.3) 

J — oo 

where ^{t — t') is determined by the correlation function (R(t)R(t')) (hereafter, angular 
brackets denote average values) [7]. Unlike the generalization Eq. (1.2), Eq. (1.3) is linear 
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in the stochastic variable x(t). Nonetheless, it is usually still difficult to solve. Often, it is 
simplified by assuming the friction to be negligible (this cannot be done if detailed balance 
is required. However, in many systems, this requirement is not necessary because the force 
is an external influence). One then ends up with the simple equation 

x = R(t) . (1.4) 

This approximation has been used recently to describe some reaction-diffusion processes 
in a force-dominated regime [8]. 

Even the simple equation (1.4) is still hard to solve when the random force is more 
complex than Gaussian white noise. In this context, solving the equation means finding the 
joint probability distribution function (pdf) for the position and the velocity, p(x, v, t), the 
marginal distribution for the velocity, p(v, t) = f p(x, v, t)dx, and the marginal distribution 
for the position, p(x, t) = J p(x, u, t)dv. Eq. (1.4) does not relate to these functions directly, 
however, and the usual approach is therefore to obtain a Fokker-Planck (FP) equation for 
the pdf's. Going from a stochastic equation like (1.4) to a FP equation can be quite 
difficult. After this is done, one still has to solve the resulting FP equation, which is far 
from trivial [9,10,11]. 

In the present paper, I outline a general approach to linear stochastic equations (in- 
cluding Eq. (1.3) and (1.4)), that often yields a closed form expression for the various 
pdf's, and may also simplify the derivation of Fokker-Planck equations. 

This approach is based on a discretization of the stochastic process represented by 
the Langevin equation by breaking up the continuous time variable into a finite number 
of steps, n. The discrete process thus generated is completely described by an n-order 
distribution function W n . In section II, I derive the fundamental formula of this paper, 
which yields p(x,v,t) as the limit n — > oo of an expression involving W n . I do not attempt 
a rigorous mathematical proof that the limit exists in general. This would be an extremely 
difficult problem since the expression obtained in section II is basically a path integral. 
However, the path integral formalism is entirely bypassed in this work in favor of a direct 
calculation approach. Indeed, the remaining sections of the paper are devoted to specific 
examples where the calculations can be carried through without using any path-integral 
method. 

Thus, in section III, the fundamental formula is applied to a class of Markov noises 
termed monovariant. This includes the Ornstein-Uhlenbeck noise, as well as, e. g., the 
Wiener and Cauchy processes. In section IV, the same formula yields a closed expression 
for the various pdf's when the noise is any Gaussian process (not necessarily Markovian) . 
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Finally, I use the general formula to derive FP-like equations: Section V treats Gaussian 
forces, and section VI treats the Markovian dichotomous noise, with which many investi- 
gations were concerned in recent years. Section VII summarizes the approach and its main 
results. 

II. General Formalism 

Eqs. (1.3) and (1.4) are particular cases of the general equation 

L(x,t) = R(t), (2.1) 

where L(x,t) is an operator acting on the stochastic variable x. I assume that L(x,t) is 
such that the general solution of the equation can be written formally as 

x = u(t)+ I g(t,t')R(t')dt' . (2.2) 
Jo 

u(t) is the solution to the homogeneous equation, and g(t : t') is the Green's function of 
the equation. Often, g(t, t') depends only on the combination t — t' ', but this is not yet 
assumed. 

Clearly, a similar relation holds for the velocity as well. Thus, we can find functions 
w(t) and h(t,t') such that 

v(t) = w(t)+ [ h(t,t')R(t')dt' , (2.3) 
Jo 

Eqs. (2.2) and (2.3), rather than Eq. (2.1), are the starting point of our approach. 
Thus, in the remainder of this paper, I assume that u(i), w(t) 7 h(t, t') and g(t, t') are given 
functions (though they may be hard to determine in practice) . For simplicity, I also assume 
throughout the paper that {R(t)} = 0. The modifications required to accommodate a non- 
zero mean value are straightforward but tedious. 

The number of stochastic variables of interest is usually determined by the order of 
the stochastic equation. Thus, in a second order process such as described by 'x = R(t), 
we have two variables of interest, x(t) and v(t). However, the method described here is 
general and applicable to any number of stochastic variables and therefore to equations of 
arbitrary order. However, for definiteness, I will describe it in the context of second order 
processes, which are the most common in physics. 

The main idea presented here is to transform the continuous process described by 
Eqs. (2.2)-(2.3) into a discrete process involving finite time steps At. Such discrete pro- 
cesses can be described easily through the n-order distribution function of the random noise 
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R(t) . It turns out that the Fourier transform of this distribution function relates very sim- 
ply to the Fourier transform of the pdf of the stochastic variable. This basic relation can 
be used to calculate the pdf explicitly in some cases, or to easily derive Fokker-Planck-like 
equations in other cases. 

To derive this fundamental relation, note that the discretized version of the random 
noise R(t) is fully determined by the set of all its n-order distribution functions 

W n (yi,ti; 2/2, ^2; - - - ; Vn,tn)dyi - ■ - dy n = Prob {y t < R(U) < y € + dy { for all i = l---n}. 

(2.4) 

i. e., W n (yi,t\; • • •; yn,t n ) dy\ • ■ ■ dy n is the probability that at the times U, the 

value of the random noise R(U) is in the range (j/i, yi + dyi). The functions W n define the 
full process R(t) through discrete and finite sets of times {U}™ =1 . This suggests describing 
the random variables x(t),v(t) through similar finite sets of times. To this end, define a 
time step At = t/n, where n is an integer that will eventually go to infinity. 

Consider now a particular realization of the process, i. e., a specific set of force val- 
ues {yi = R(iAt)}2 =1 . For this realization, define sets {£i}™ =1 and {i>i}?=i through a 
discretization of Eqs. (2.2) and (2.3), i. e.,: 

n 

x n = u n + y^ j g(n, j) yjAt , (2.5a) 

n 

v n = w n + ^2h(nJ)yjAt, (2.56) 

i=i 

where 

u n = u(nAt) , w n = w(nAt) , (2.5c) 

g(nj) = g(nAt,jAt) , h(n,j) = h(nAtjAt) . (2.5d) 

As n — > 00 and At — > with nAt remaining constant, the sets {xi} and {vi} converge to 
the continuous processes x(t), v(t). 

Let us now define a discrete analog of the joint pdf p(x,v,t). This discrete pdf, 
denoted p(x,v,n), is defined as 

p(x, v, n) dx dv = Prob {x < x n < x + dx and v < v n < v + dv} , (2.6) 

where x n and v n are now defined by Eqs. (2.4) and (2.6), with respect to the values 
of the random noise R(t). Our aim is to calculate p(x,v,n) explicitly and hope that as 
n — > 00, At — > 0, p(x, v, n) will tend to p(x, v, t). 
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The probability p(x, v, n) is the total probability of all possible realizations {yi} such 
that x < x n < x + dx and v < v n < v + dv, where x n and v n are given by Eqs. (2.4)-(2.6) 
and n = t/At. For simplicity, I assume that the set {yi}™ =1 may take any value in the 
range (—00, 00) (if the range is discrete, one can use appropriate 5-functions to restrict the 
values of yi). Therefore, we have 

/oo roo 
dyi--- dy n W n (y 1 ,At; y 2 ,2At; . . .; y n ,nAt)S(x-x n )S(v-v n ) . (2.7) 
-00 J — 00 

The two 5-functions in the integral, S(x — x n ) and S(v — v n ), ensure that only processes 
{Ui}7=i which yield the appropriate final positions and velocities are counted. Note that 
x n and y n are functions of {yi} through the relations (2.5). For simplicity, I shall write 
from now on W n (y 1 , ...,y n ) instead of W n {y u At; y 2 , 2At; ...;y n , nAt). 
Eq. (2.7) can be simplified by introducing the Fourier transforms (FT), 

/oo />0O 
dx / dve i{xk+vd) p(x,v,n). (2.8a) 
-00 J — 00 

/oo />oo 
dyi--- dy n exp [i(yi<p! H h y n <Pn)] W n (y 1 , ...,y n ) (2.86) 
-00 J — 00 

Substituting Eq. (2.7) into Eq. (2.8) and interchanging the order of integration, we 
obtain: 

/oo />oo />oo />oo 

dyi--- dy n dx dv e [xk+v0) W n {y u . . . , y n ) 6(x - x n )S(v - v n ) . 
-00 J— 00 J— 00 J— 00 

(2.9) 

Performing the integration over x and v and replacing x n and f n with the corresponding 
expressions from Eqs. (2.5) which relate them to the values y i: we have: 



/oo poo 
d Vl --- dy n e^ k+w ^W n (y u ...,y n ) 
-00 J —00 



x exp jz^[%(n,z) + Oh(n, i)]yiAt j . (2.10) 

Defining <pi = [o(n, i)k + h(n, i)0]At, and comparing with Eq. (2.8b), we see that 

p(k, 9, n) =e l ( Unk+w ^W(vu ...,<p n ), (2.11a) 

where 

L p i =[g(n, i)k + h(n, i)9]At . (2.116) 
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Ultimately, however, we are interested in the FT of p(x, v,t), 

/oo />oo 
dx / dve^ xk+v0 ^p(x,v,t). (2.12) 
-oo J — oo 

To obtain this, we take the limit n — > oo, At — > of 6*, n), i e., 

p(fc, 61, t) = Jim e *-*-H»»«) wfyx, ...,?„), (2.13a) 

At^O 
nAt—t 

with 

(p t = [g(n,i)k + h(n,i)6]At. (2.136) 

Eq. (2.13) is the central result of this section and the basis for the rest of this paper, 
which deals with specific applications. 

Because of its importance, one should note that while reasonable from the physical 
point of view, Eq. (2.13) nonetheless raises some subtle mathematical questions. We 
cannot prove that the limit n — > oo, At — > of p(k, 8,n) exists in general, nor that given 
its existence, it is indeed equal to p(k,6,t). The source of this difficulty is Eq. (2.7), 
which is clearly a discrete path integral, and therefore subject to various mathematical 
reservations when the limit n — > oo is taken formally. 

However, in this paper the usual formalism associated with path integrals is bypassed 
in favor of the direct definition Eq. (2.13). In other words, we do not take the limit 
n — > oo, At — > formally and write the resulting expression as a path integral. Rather, 
all calculations are carried out for the finite n case, and the limit then taken directly. It 
is remarkable that in many cases this limit turns out to be well defined mathematically 
as well as calculable, with the only assumptions being the smoothness properties usually 
postulated in physical problems. The next two sections are devoted to such calculations 
and cover a wide range of cases. 

III. Monovariant Markov Noise Processes 

As a first example where the limit in Eq. (2.13) can be calculated explicitly, consider 
a class of Markov noise processes defined below: 

A Markov process is fully determined by the initial distribution P (yi) : and the tran- 
sition probability, T(j/j + i,tj + At\yi,ti), from the value yi of the random function R(U) 
at time U to the value yi+\ = R{U + At) at a time U + At [12]. It is usually assumed 
that this function depends only on the time step At and not on tj. I shall assume this for 
simplicity, but the calculations can be modified for more general cases. 
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I shall say a Markov process is monovariant if the transition probability does not 
depend on y i+ i and x/i separately, but rather on a single linear combination of these 
variables. Such a combination can always be written as yi + i — //(At)j/j, where //(At) 
is some function of the time step At, and I also assume that //(0) 7^ 0. Hence, 

T(y i+1 ,At\yi) = T At (y i+1 - n(At)y t ) . (3.1) 

For any Markov process, the functions W n (yi, . . . , y n ) are given by [12]: 

Wniyx,..., y n ) = P (yi) T At (y 2 \yi) T At (y 3 \y 2 ) ■ ■ -T Af (y n |// n _i) , (3.2) 

where we assume all transitions take place during a time-step At = t/n, with n eventually 
going to infinity. 

Many of the noises used in physical problems are monovariant. Some examples are the 
Wiener noise (with //(At) = 1), the Cauchy noise (also with //(At) = 1), and the Ornstein- 
Uhlenbeck (O-U) noise, for which //(At) = exp(— At/r), where r is the correlation time 
[12]. The O-U noise is widely used to describe forces with finite correlation time (e. g., 
Ref. [9]). The r — > limit corresponds to white noise. 

Let us now calculate the pdf for a process driven by a general monovariant Markov 
noise: The multi-variable FT of W n {y\, . . . , y n ) is 

/oo />oo 
dyi ■ ■ ■ dy n exp [i(yi<pi H h y n <fn)} Po(yi) 
-00 J — 00 

x T At (y 2 \ yi ) T At (y 3 \y 2 ) ■ ■ ■ T At (y n |y n _!) . (3.3) 

If the Markov noise is monovariant, so that T At {y'\y) = T At (y' — //(At)//), we can perform 
the integration in the following way: Define a single- variable Fourier transform T such 
that 

/oo 
e*" T At ( V ) dn , (3.4) 
-00 

where rj = y' — //(At)// is the single variable on which T A t(y'\y) depends, except for any 
direct dependence on At. For conciseness, I assume such a dependence implicitly and do 
not write the At subscript anymore. We now change the integration variables {yi}f=i in 
Eq. (3.3) to the new set: 

m=yi> ( 3 - 5 ) 

V2=V2~ n{At)y x , 
Vn=y n - //(At)// n _i . 



8 



The Jacobian of this transformation is clearly unity, and each T^t (yi\yi-i) depends on rji 
only. To rewrite Yl7=i VWi m terms of the new variables, we invert Eq. (3.5), which gives: 

Vi=Vi, (3-6) 
V2 = V2 + MA£)t7i , 

V3 = vs + MAt)??2 + [MAt)] 2 m > 

i 



Hence, 



k=i 



n n ( n 

1 k — i 



Y VWi = Y ViS Y W**)\ 1 \ . (3.7) 



i=l i=l (. k=i 



Substituting Eq. (3.7) into the expression for W, Eq. (3.30), and using the definition of T, 
Eq. (3.4), we finally obtain 



W^ l ,...^ n ) = P [C]l[f(C l ) , (3.8a) 



1=2 



with 



C = J> fe [MA*)] fe_1 , (3-86) 
k=i 

n 

& = ma*)]*-* ■ (3.8c) 



k=i 



In Eq. (3.8a), 



/oo 
e« v P (y)dy. (3.9) 
-oo 



In accordance with the general formula (2.11), p(k,9,n), is obtained by replacing <pi 
in Eqs. (3.11) with [g(n, i)k + h(n, i)0]At. This yields: 

n 

p(k, 9, n) = e^ k+w ^P [C] J] f (&) , (3.10a) 

i=2 

where 

/ n 

iife-1 Ail , /i/V^l/ nr ui\iH 



C = ^^^^)[MA*)] A^J +^^/i(n,fc)[MAt)] fc - 1 Atj , (3.106) 
E </(n, fc) HAt)] k ~ l At)+e[J2 Kn, k) HAt)] k ~ l At) . (3.10c) 



^k=i / \k=i 
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Last, we must take the limit n — > oo, At — > 0. This is a somewhat tedious step, and its 
details can be found in Appendix A. The final result is the following very general formula: 



where 



B(0 = 



d 



dAt 



logT At (0 



J At=0 



C(t) = k I g(t,s)V(s)ds + I h(t,s)W(s)ds , 
Jo Jo 

^(t,t') = k I g(t,s)V(s-t')ds + [ h(t,s)V(s -t')ds , 
Jt' Jt' 



with 



^(z) = k(0) exp[a^ 
d[\ogv(t)} 



a 



dt 



1 d/i(t) 



(3.11a) 

(3.116) 
(3.11c) 

(3.11rf) 

(3.11e) 
(3.11/) 



t=o 



Eqs. (3.11) are general and can be applied for a great variety of noise processes and 
of dynamics (i. e., forms of the Langevin equation). The price of this flexibility is the 
complex form of the equations. Their meaning and use can be clarified, however, by using 
as an example a very simple equation which was investigated recently [8,9], 

x = v = R(t) , (3.12) 

The solution is v = vo + J * R(t')dt', so that we have (compare Eq. (2.5)): 

w(t) = vq , (3.13a) 
h(t,t') = l. (3.136) 

t t t' 

Now x(t) = xq + J Q dt" v(t") = xq + vot + J Q dt' J Q dt" R(t"). According to a well-known 

theorem [13], 

[ dt' [ dt" R(t") = [ dt' (t - t')R{t') , 
Jo Jo Jo 

so that 



(3.14) 



u(t) = x + v t , 
g(t,t')=t-t' (t'<t). 



(3.15) 
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Let us now take R(t) to be an O-U noise, defined by [12, 9]: 

,2 



Poiyi) = 
T&t{yi+i\yi) = 



T 



7r/ 2 [l-^(At)] 



At 



exp 



/o 2 [l-^ 2 (At)] 



//(At) = exp . 



where fo is a parameter describing the intensity of the random force. 
The Fourier transforms P and T, are (see Eqs. (3.4) and (3.9)) 

i 5 o[C]=exp (- f le 



At 



TAt(0 = ex P 



(l-^ 2 )/o 2 
4r 



e 



From Eq. (3.1 If), we find 



a = — , 

T 



so that is now (from Eq. (3. lie)) 

V(z) = exp 

From the definition of B(£), Eq. (3.11b), we have: 



where 



«(*) = / 2 

m = \fo 
i{t) = ft 



j - U 2 r + r 3 (l - e"^) - rhe-^ 

-rt(l- e~ th 
t-r(l- e~ t/T 



(3.16a) 
(3.166) 
(3.16c) 



(3.17a) 
(3.176) 

(3.18) 

(3.19) 

(3.20) 



The details of the remaining algebra are in Appendix A. The final result is 

p(k, 9, t) = exp { % [k(x + v t) + 9v ] - - [a(t)k 2 + 2(3(t)k6 + -f(t)6 2 ] 1 , (3.21a) 



(3.216) 

(3.21c) 
(3.21d) 



Since this is a Gaussian Fourier transform, it can be inverted easily to yield the real 
space joint pdf: 

1 

p{x,v,t) = 



2ttV(«7-/9 2 ) 



x exp 



2(~Z b( x ~ x o~ vot) 2 - 2(3(x -xq- v t)(v - v ) + a(v - v ) 2 ] J , 

(3.22) 
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with a, (3 and 7 as denned in Eq. (3.21). 

The marginal distributions for the position, p(x, t), and for the velocity, p(v,t), are 
most easily obtained from Eq. (3.21) by setting k = and 9 = respectively, then Fourier 
inverting. The results are: 



p(x,t) = 
p(v,t) = 



1 



y/2ica(t) 
1 



exp 



exp 



(x — Xq — Vot)' 
2a(t) 

27(f) 



(3.22a) 
(3.226) 



Eq. (3.11) with an O-U noise was investigated by Heinrichs [9], who derived Fokker- 
Planck (FP) equations for the joint pdf as well as for the two marginal pdf 's p(x, t) and 
p(v, t). He solved these equations approximately only in two limiting cases: For very small 
t (near r = 0, which corresponds to white noise) and for very long r (near 1/r = 0, which 
correspond to a constant force). 

Here, however, we have obtained exact results in Eqs. (3.22) and (3.23), quite easily 
and without solving any differential equation. These expressions reduce to Heinrichs's 
results in the appropriate limits. Also, a simple substitution shows that Eqs. (3.22)-(3.23) 
are indeed solutions of Heinrichs's FP equations [9] . This further confirms the effectiveness 
of the present approach. 



IV. Gaussian Processes 

Another class of random noises that can be solved completely is the group of Gaussian 
noises. Except for the O-U noise, all these processes are non-Markovian, and they are 
defined by the requirement that all the distribution functions W n (yi, . . . , y n ) be of Gaussian 
form. For zero- mean noises, this is equivalent (e. g., Ref. [7]) to the requirement that their 
Fourier transform be given by: 



W((pi,...,(p n ) =exp 



1 n 



(4.1) 



where 4>{i,j) is the correlation function of the noise at time steps i and j, i. e., (yiVj). 
Note that formula (4.1) assumes that 4>(i,j) = 4>(j,i)- 
Substituting (4.1) into (2.11), we have: 



p(k, 9, n) = exp<( i(ku n + 9w n ) - ^ 



k 2 Yl 9(n,i)g(n,j)(j>(i,j) Ar 

■ i,j=0 
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+k9 [g(nJ)h(nJ) + h(n,i)g(nJ)](P(iJ) At 2 

(4.2) 



i,j=0 
n 



+6 2 Hn,i)h(n,j)<j)(i,j) At 

i,3=0 

Upon taking the limit n — > oo, At — > 0, all the sums become integrals. Hence, we have 
from Eq. (2.13): 

p(k,e,t) = exp\i(ku(t) + 6w(t))- - [a(t)k 2 + 2(3(t)k6 + -f(t)6 2 ]\ , (4.3a) 



where 



a(t)= f dt' f dt" g(t : t')g(t : t")(j)(t\t") : (4.36) 

= 1 /V/V [g(t,t')h(t,t") + h(t,t')g(t,t")} <t>(t',0 , (4.3c) 
^ Jo Jo 

/ dt' f dt" h(t, t')h(t, t") (f)(t', t") , (4.3d) 
Jo Jo 



7(t) 



and <j)(t',t") = (R(t')R(t")) is the continuous time correlation function of the noise R(t). 

As an application of Eq. (4.3) and as a consistency check, we can use Eq. (3.12) again, 
x = R(t), with R(t) as O-U noise (since this noise is both Markovian and Gaussian). For 
this case 



2 



W -t"\ 



(4.4a) 

g(t,t')=t-t' (t>f), (4.46) 

h(t,t') = l. (4.4c) 

The functions a, j3 and 7 calculated through Eq. (4.3) are identical (as expected) with 
those obtained for the same noise (O-U) in Eq. (3.21). 

V. Differential Equations for PDF's in the Gaussian Case 

The preceding sections presented cases where direct calculation of the pdf could be 
performed, i. e., when the limit in Eq. (2.13) could be calculated explicitly. In some cases, 
however, the limit proves too complicated for explicit calculation. Eq. (2.11) then offers an 
alternative point of view by providing a general approach for deriving FP-like differential 
equations for the pdf's. 

To exemplify this approach through a simple example, consider the Gaussian random 
noises described in section IV. Having arrived at Eq. (4.2), we now choose not to take the 



13 



limit n — > oo,At — > which would yield Eqs. (4.3). Instead, let us calculate the time 
derivative of p(k, 6>, t) from its discrete form p{k, 6>, n). This is just 



dp(k,9,t) ,. p(k,9,n+ 1) — 0, n) ,. Ap 

= hm = hm — — 

dt At^o At At^o At 



(5.1) 



With the help of Eq. (4.2) for p(k, 9, n), this calculation is straightforward but tedious, 
and the details can be found in Appendix B. The final result is: 



dp(k,9,t) f.d-u, .dw n 1 

k = Si—rk + i— —9 

dt I dt dt 2 



dt 



da 
dt 



dt 



p(k,9,t), (5.2) 



where ck,/3 and 7 are the functions defined in Eqs. (4.3). 
Eq. (5.2) can be Fourier inverted to yield: 



dp(x, v, t) du dp dw dp 1 
dt dt dx dt dv 2 



da d 2 p ^ d(3 d 2 p d'j d 2 p 



dt dx 2 



dt dxdv dt dv 2 



(5.3) 



Here again, one should note the high generality of these equations. The result obtained 
here is valid for a wide range of Langevin equations as well as for a wide range of random 
noises. 

We can now apply this formalism to the simple equation x = R(t) with the random 
force being the O-U noise already considered in section III for which we have found the 
pdf [Eq. (3.22)]. Then from (4.3) and (4.4), we have 



~dt 

d£ 

~dt 

da _ 2 

~dt~ h 



1 

2 



1 - 

f 2 \t(2 - e -*/r) - T (l - 



t 2 - tr(l 



-t/n 



(5.4a) 
(5.46) 
(5.4c) 



The remarkable thing about Eqs. (5.3)-(5.4) is that they differ from the FP equa- 
tions derived by Heinrichs for the same stochastic process [9]. Instead of Eqs. (5.3)-(5.4), 
Heinrichs obtains (for the case he investigated, i. e., xq = vq = 0, hence u = w = 0) : 



dp(x, v, t) 



dt 



d 



d' 2 



V dx dxdv 



+ a(t) 



8P_ 

dv 2 



p{x,v,t) , 



where 



-t/r 



ait) = f 
b(t) = £[(t + T) e"'/' 



(5.5a) 

(5.56) 
(5.5c) 
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This equation is not wrong. Indeed, as mentioned before, the function p(x,v,t) calcu- 
lated in Eq. (3.29) is the solution of Heinrichs's equation (5.5), as well as the solution of 
Eq. (5.3) (both with the appropriate initial conditions). This means that the FP equa- 
tion corresponding to a given Langevin equation is not unique (and the existence of two 
different equations means trivially that there is an infinite number of them). Rather, the 
form of the equation depends on the specific method used to derive it. All the various FP 
equations, however, are equivalent to each other through the specific mathematical form 
of the solution. 

It may seem that this last observation, as well as having the explicit solution, renders 
all discussions about the equations unimportant. This may be true for the case of the 
unconstrained particle, which is discussed here, but the discussion is very relevant to cases 
of constrained particles. 

A particle can be "constrained", e. g. : by the addition of boundaries to the system (for 
example, absorbing boundaries are often added for calculating first passage times). In the 
theory of the classical FP equation [14], the presence of such boundaries doesn't change 
the equation. Rather, it implies new mathematical boundary conditions onp(x, v, t) (these 
usually make the task of solving the equations much more difficult). 

It is tempting to extend such reasoning to Fokker-Planck-Zz/ce equations, such as 
Eqs. (5.3) or (5.5). Indeed, Heinrichs used precisely this approach to tackle the case of a 
particle driven by O-U noise when two absorbing boundaries are present [15]. However, 
the non-uniqueness of the FP-equation casts serious doubt on the correctness of such an 
approach. It is highly unlikely, for example, that Eq. (5.3) and (5.5) will still be equivalent 
to each other when new boundary conditions are added to them. Furthermore, it is quite 
possible that neither of them is correct in the presence of boundaries. Indeed, Eq. (5.3) 
follows from the general Eq. (2.11), which assumes implicitly the absence of constraints in 
so far as it allows summation over all paths that lead to the adequate final position and 
velocity. On the other hand, Heinrichs derived Eq. (5.5) by applying some averaging iden- 
tities to the stochastic differential equation (3.12) [9]. Under constraints, such averaging 
procedures may change and one cannot assume a priori that the FP equation remains the 
same (or indeed exists at all). 

In other words, the non-uniqueness of the FP-equations implies that whenever we add 
boundaries to the problem, we must provide independent mathematical justification for 
any FP equation we choose to use. In the absence of such justification, the results obtained 
by adding boundary conditions to any FP equations cannot be trusted. 

VI. Differential Equations for PDF's in the Case of Markovian Dichotomous Noise 
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For Gaussian noises, we have both an explicit expression for the pdf and a FP-like 
equation Eq. (5.3). However, and even more importantly, the method described in the last 
section can be used in cases where the limit n — > oo, At — > proves too complex to perform 
explicitly. As an example of the usefulness of this approach, I discuss now a specific type 
of random driving force, the Markovian dichotomous noise, which is used, for example, in 
the theory of second order Butterworth filters [16]. 

Masoliver [10,11] has treated two stochastic differential equations with this particular 
noise, i. e., x = R(t) and x + "fx = R(t), and has obtained FP-like equations for p(x, v, t), 
p(x, t) and p(v, t) in both cases. Deriving the equations requires much work in each case, yet 
the general formula (2.13) provides a unifying point of view from which all these equations 
follow naturally. In this section, I derive a generic differential equation for all the various 
pdf's, when the driving force is the Markovian dichotomous noise, but the form of the 
equation is still general. Masoliver's equations follow immediately as several particular 
cases of this general equation. Thus, we shall obtain a powerful generalization as well as 
a compact form for all the cases studied. 

The Markovian dichotomous noise is defined by: 

PoM = \ [%i - a) + 6( V1 + a)] , (6.1a) 

TaMv') = \ [f(At)S(y - y') + g(At)S(y + y')} , (6.16) 

where 

f(At) = 1 + exp (-2AA£) , (6.1c) 

g(At) = 1 - exp (-2XAt) . (6.1d) 

Thus the force alternates between two values only, a and —a, and the probability that a 
switch from one value to the other occurs in the interval (t, t + dt) is \e~ xt dt [10]. The 
parameter A is the average time between switches. 

The calculation of p(k : 9, n) follows the lines described in sections II-IV. The details 
are in appendix C. The final result is 

p(k, 9, n) = exp [i(ku n + 9w n )] E ( 2 ) / n_fe_ V _1 

x [F^i, ...,<p n ) + FZi-n, -ip n )] , (6.2a) 

where 

(Pi = [kg(n,i) + 9h(n,i)]At. (6.26) 
The coefficients Fj}((fi, . . . , <p n ) are determined recursively by 

F2(<p 1 ,...,<p n ) = e ia * 1 [^K.-.^ + d 1 ^,.-.,-^)] , (6.3a) 
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with the conditions 



F™(<p u ...,<p n ) 



0. 



0. 



(6.36) 
(6.3c) 
(6.3d) 



In principle, Eq. (6.3) is the full solution of the problem. However, the limit n — > oo, At — > 
cannot be taken explicitly in all generality. Nonetheless, Eq. (6.3) provides an approxima- 
tion of p(k, 6>, £), which can be calculated by selecting a sufficiently high value of n. 

In the following, however, I use Eq. (6.3) for deriving the differential equation obeyed 
by p{k, 6>, t) along the lines of section V. In the present case, however, this procedure yields 
a second order equation in the time derivatives, rather than a first order one as in the case 
of Gaussian noise. 

For the purpose of this derivation, I shall assume that g(n,i) and h(n,i) depend on 
the difference n — i only (or, in the continuum version g(t, t') = g(t — t')). This is usually 
true. Also, I assume that u{t) = w(t) = (the particle starts at rest from the origin). 
This, however, is no loss of generality, because if we are interested in other cases, all we 
need to do is multiply the function p(k,9,t) obtained for the case u(t) = w(t) = by 
exp[i(ku(t) +6w(t)]. 

The details of the calculations can be found in appendix D. The main results are: 
There is no first order equation in the time derivatives. This is because, to fist order 
in At, we find that 



Ap(k, 6, n) = ion £ ( - J /»-*- V" 1 [F£( V 2, • • • , ¥>n+i) - F£(-<p 2 , . . . , -<p n+1 )} , 



with cpi having the same meaning as in Eq. (6.2), i. e., tpi = [kg(n, i) + 6>/i(n, i)] At. Com- 
paring this expression with Eq. (6.2) reveals that Ap cannot be expressed in terms of 
p(k, 9, n) only, hence we cannot write down a first order equation. 

We therefore look at the second derivative, for which we need to calculate A 2 p = 
p(k,9,n + 2) + p(k,9,n) — 2p(k,9,n + 1). This expression is then expanded to second 
order in powers of At. As shown in appendix D, we can then express A 2 p/At 2 in terms 
of Ap/ At and p, so that in the limit n — > oo, At — > we obtain a second order differential 
equation for p(k, 6,n). Quoting from appendix D, we have finally that 



n 




k=l 



(6.4) 



d 2 



] p(k,9,t) 
dt 2 




1 dif\ dp(k,6,t) 
ip dt J dt 



+ a 2 ip 2 p(k,9,t) = 0, 



(6.5a) 



where 
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(p = kg(t) + 9h(t) . (6.56) 

The power of Eq. (6.5) lies in the generality of the functions g(t) and h(t) that appear 
in it. This is because Eq. (6.5) covers all linear stochastic processes driven by a Markovian 
dichotomous noise. Thus, although the noise is fully determined, the dynamics remains 
arbitrary [within the confines set by Eq. (2.2)], and therefore, Eq. (6.5) covers at once a 
wide range of stochastic processes. 

Consider, for example, the simple equation x = R(t), for which g(t) = t,h(t) = 1. 
Then , from Eq. (6.5) we have immediately for p(k,9,t): 

d 2 p%e,t) f _ _fc_\ dp(k,e ,t) , ^ ^ , m2 

dt 2 V kt + ej dt 

The equations for the FT of the marginal distribution for the position, i. e., p(k,t) and 
for the velocity, i. e., p(9,t) follow at once by setting 9 = and k = respectively: 

92p{k > t) + (2\-l) + aW m t) = o . ( 6 .7a) 



+ " 77^ ! + « 2 (** + ef p(k, 0,t) = 0. (6.6) 



dt 2 \ t dt 



Eq. (6.7b) is the FT of the telegraph equation for p(v,t). The two Eqs. (6.7) were derived 
separately by Masoliver [10], along with a somewhat more complicated version of Eq. (6.7) 
(the complication being due mainly to notations). As Masoliver pointed out, k and 6 
do not appear independently in the equation for p(k,9,t), but do so only through the 
combination kt + 9. This fact was not immediately apparent in his version of the equation 
for p(k, 6,t), but it was used later to derive Eq. (6.6). 

There is no obvious reason in Masoliver's work for the appearance of such a combina- 
tion. Here, on the other hand, it appears as a natural consequence of the general formula, 
Eq. (2.13), in which k and 9 are always bound in a well defined combination (except where 
they relate to u(t) and w(t)). This also makes very clear the meaning of the particular 
coefficients in this combination (t and 1), as they are merely the Green's functions of the 
original stochastic equation. 

Masoliver later tackled the more general equation x + = R(t) [11]. Much work is 
needed to derive the FP equations for p(k, 9, £), p(k, t) and p(9, t). Moreover, in Masoliver's 
approach, the equations for the marginal distributions do not follow immediately from the 
equation for p(k, 9, t), but instead require still more work. In the present approach, on the 
other hand, the derivations are almost trivial. For the equation x + "fx = R(t), we have 

g(t -t') = -{l- exp [- 7 (t - £')]} • (6.8a) 

7 

h(t-t') = exp[- 7 (£-£')] • (6-86) 
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(easily derived, e. g., from Laplace transforming the equation). Thus, we have at once: 



d 2 p(k,9,t) 
dt 2 



+ 



2A 



7(fc-7^)e~ 7t 
k - (k - 7#)e" 



■7* 



dp(k,9,t) 



+ ^ [k - (k - 7 6»)e- 7t ] z p(fc, 6»,t) = . (6.9a) 



|2 



d 2 p(M) 

at 2 



+ 



2A 



7e" 



- 7 t 



1 - e-T* 



dp(k,t) a 2 k 2 



+ — (l-e" 7t ) p(k,t) = 0. (6.96) 



d 2 p(e,t) 

dt 2 



+ (2A + 7 ) 



at 7 
dp(e,t) 



dt 



+ a 2 e 2 e~ 2 ' rt p(8,t) = 0. (6.9c) 



These agree with Masoliver's equations, but are derived much more easily. Again, the fact 
that k and 9 appear only through the combination k [1 — exp(— jt)] + O'j exp(-'jt) is not 
evident in Masoliver's approach, while it finds a very natural explanation in the present 
framework. 



VII. Summary 

The main thrust of this work is methodological and conceptual. Rather than con- 
centrating on a specific process, I have presented a general approach to a large class of 
stochastic differential equations. 

The basis of the method is the discretization of the process described by the equation 
into n time-steps. Each realization of the process consists then of n steps, and its probabil- 
ity can be calculated from the elementary properties of the random noise. In this discrete 
approximation, the probability of finding the particle in a given state at time t is the sum 
of the probabilities of the realizations which lead to this final state. Taking n — > oo then 
yields the required probability for the continuous process. 

This method is successful on at least two fronts. First, if the limit n — > oo can be 
calculated explicitly, we obtain a closed expression for the required probability distribution 
function. This turns out to be the case for all Gaussian processes and monovariant Markov 
processes. There is no reason to believe these exhaust the possibilities of this method. 
Quite likely, other processes can be completely solved in this way. 

Second, the expressions obtained for the discrete approximation are useful for de- 
riving Fokker-Planck equations for the probability distribution functions. For the case 
of the Markovian dichotomous noise, I have shown that there is a generic second order 
FP equation, which covers all the possible stochastic differential equations to which the 
general method is applicable. Thus, we obtain a unified point of view on the various FP 



19 



equations this particular noise can generate. In particular, for a joint probability distribu- 
tion involving more than one stochastic variable (e. g., position and velocity), the Fourier 
frequencies corresponding to the stochastic variables appear only in a well specified linear 
combination. The coefficients in this combination are the Green's function attached to 
each stochastic variable. Thus, an important mathematical property of the FP equation 
is made explicit and its origin is clarified. This was not obvious in other approaches which 
have been applied to this problem. 

Another general point underscored by the present approach is that there may be 
several FP equations corresponding to a single stochastic differential equation. This non- 
uniqueness may have important implications when we consider particles constrained, e. g., 
by various boundaries. The standard recipe in Fokker-Planck theory, i. e., adding some 
boundary conditions, is incomplete, since the various FP equations may yield different 
results when supplemented with the same boundary conditions. Thus, we must be careful 
when applying FP equations derived for free particles to constrained cases. Justification 
is always needed when taking such a step. 

Finally, the discrete approximation should be of interest for numerical estimates of the 
pdf. I have not tried to address this issue here, and concentrated on analytical applications 
instead, but this should be investigated by interested parties. 

In the last few years, there seem to be somewhat more interest in stochastic equations 
with complex environments and simple noises than in equations with simple environments 
and complex noises. I feel at least part of the reason for this is that the latter cases have 
proved quite difficult mathematically. Certainly, many physical systems do exhibit complex 
random noises. The method presented here may simplify the mathematical difficulties to a 
considerable extent. I hope this will motivate interested researchers to use it and develop 
it further. 
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Appendix A 

We wish to take the limit n — > oo, At — > of Eq. (3.10), i. e., 

n 

p(k, 9, n) = e^+^Po [C] J] f (6) , (Ala) 

where 

c = k g(n, k) um k ~ l *tj + o (j2 h ^ k ) H^)\ k ~ l A ^ > ( A - lh ) 

^i = k(j2 k ) H^t)] k ~ l AH + 9 [J2 Hn, k) [^(At)] k - 1 At) . (A.lc) 
\k=i J \k=i J 

To make this step clearer, I'll perform it in several stages. 

First, the expression [n(At)] k ~ % seems to become ill-defined in this limit. Let us 



therefore introduce the following quantity: 

„_d[log //(*)] 



dt 



1 dn(t) 



,o Kt) dt 



(A.2) 



Hence, to first order in At, 



fi(At) « (i(0) exp(aAt) , (A3) 

so that 

[/j(At)] k ~ l ~ ^(0) exp [a(k - i)At] = (Jfe - i) . (A4) 

This expression has a well defined value in the limit n — > oo, At — > 0, when the function 
\1/ becomes 

= //(0) exp(az) . (A5) 

Next, note that all sums of the form X)fc=i ' ' ' At will become integrals in this limit. 
Introducing an integration variable s = kAt, an initial time t' = zAt and the final time 
t = nAt, and using the definition of Eq. (A. 5), we have: 

V <,(n, fc) ^(At)] fc " 1 At — / git, s)V(s- t') ds . (A6) 



k=i 



Similarly, 



J2h(n,k) [n(At)} k ~ l At — ► / h{t,s)y{s-t')ds, {Ala) 

V g(n, k) [^(At)]"" 1 At — / s ) *( s ) ds , (A 76) 
fc=i Jo 

V/i(n,fc) [^(At)] fc_1 At — > / s) (s) ds , (A7c) 



fc=i 
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where in the last two lines, the difference between (k — l)At and kAt is neglected, since it 
yields a vanishing correction of order At to the expressions on the r.h.s. 

Finally, consider the product over all the T-functions in Eq. (A. la). Following the 
idea used in Eqs. (A.2)-(A.5), we define 



5(0 



d 



dAt 



logf At (0 



At=0 



(A.8) 



From the definition of TAt(£)> we have Ta*(0) = 1. Hence, in the limit n — > oo, At — > 0, 



T At (0»exp[S(0At] . 



(A9) 



Substituting Eq. (A. 9) into the expression n"=2^(^) which appears in Eq. (A. la) 
yields 



i=2 



i=2 



(A10) 



In the limit n — > oo, At — > 0, we have, from Eq. (A.lc), that 



& — >^(t,t') = k [ g(t,s)y(s-t')ds + 9 ! h(t,s)y(s-t')ds, (All) 



and therefore 



where 



«■ ,>t 

i=2 ^ 

Putting together Eqs. (A.l), (A.6), (A.7) and (A.12), we finally have that 
p(k, 6, t) = e i[ku(t) + e w (t)]p Q m] exp y* dt , B m t , )} | ^ 

C(t) = k I g(t,s)V(s)ds + I h(t,s)W(s)ds, 
Jo Jo 

pt pt 

^(t,t') = k g(t,s)V(s-t')ds + 6 h(t, s)#(s - t') ds , 
Jt' Jt' 



(A12) 



(Al3a) 

(A136) 
(Al3c) 



As an example, let us now apply this general result to the specific case of x = R(t). 
As noted in Eq. (3.15), 



u(t) =x + v t , w(t) = v , 
g(t,t') =t -t' , h(t,t') = l 



(A14) 



(f < f ) 
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Let us now take R(t) to be an O-U noise, denned by [12, 9]: 
PoiVi) = 
TAtiVi+ilVi) = 



—72 ex P 



T 

72'. 



exp 



7T/ 2 [l-^(Af)] 

//(At) =exp • 

where /o is a parameter describing the intensity of the random force. 
The Fourier transforms P and T, are (see Eqs. (3.4) and (3.9)) 



Tai(0 = exp 



4r 



:i-^ 2 )/o 2 , 

~Tr ' 



From Eq. (A.2), we find 



a — — 

T 



so that ^(z) is now (from Eq. (A. 5)) 

ty(z) = exp • 
From the definition of B(£), Eq. (A. 8), we have: 



5(0 = -§e 



After some algebra, we obtain 



/ h(t, s)*(s - t') ds = I exp 
Jv Jv 

[ g(t,s)ty(s -t')ds = [ (t-s)exp 
Jv Jv 



(s-t>) 



ds = r < 1 — exp 



(t-f 



- s) exp 


- ( a -f)l 


ds 




T 






t-t 


+ r ^exp 


r (*-*y 


-)] 


r 


r 



from which we have immediately (by setting t' = 0) 

/ h(t,s)W(s)ds = T 
Jo 

[ g(t,s)*(s) 

Jo 



1 — exp — 



(is = T it — T 
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1 — exp — 



Since ( = k J Q g(t, s)*(s) ds + 9 f Q h(t, s)*(s) ds, by substituting (A.21) into (A.16a) we 
obtain 



P (C) = exp j-^ [a(t)/c 2 + 6(t)# 2 + 2c{t)k$] j , 



where 



a(t) 
b(t) 
c(t) 



ft 
2 

ft 
2 



t - T 



1 — exp ( — 



t 



1 — exp ( 



1 — exp ( — 



t 



— T 



1 — exp ( — 



(A22a) 

(A226) 
(A22c) 
(A22d) 



Similarly, we must substitute Eqs. (A. 19) and (A. 20) into the expression J* dt'B t')] , 
where £(t, *') = fe s)*(s - ds + 9 f* h(t, s)*(s - £') (is . The integral over t' can 

be performed explicitly, and after some more algebra , we have that 



/' 

Jo 



dt'B[Z(t 1 t')] = 



- [l(t)k 2 + m(t)9 2 + 2n(t)k9] , 



with 



l(t) = ft 
m{t) = f 

n(t) = ft 



1 +tT 2 -t 2 r + ^-- 2tT 2 e~ t ' T - T -e~ 2t ^ 

O Zi _ 



t--r + 2re~ t / T - -e" 2 */ 7 
2 2 



t --tr+ T -- r 2 e-^ + r 2 e- 2 */ T + tre^ 



(A23a) 

(A236) 
(A23c) 
(A23d) 



Substituting (A. 14), (A. 22) and (A.23) into (A. 13) yields the final answer: 



p(k, 9, t) = exp \ i [k(x + v t) + 9v ] - - [a(t)k 2 + 2[5(t)k9 + ~/(t)9 2 ] \ , (A24a) 



where 



<*(t) = ft 



3 2 



t 2 r + r 3 (l - e-*/ T ) - r 2 te-^ 



m = \ft\t 2 ~rt (l-e-^) 

7(t) = ft [t-^i 1 - e_t/r ) 



(A246) 

(A24c) 
(A24d) 



Appendix B 

To calculate the discrete form of the time derivative of p(k,9,t), we first expand 
p(k, 9, n+ 1) in powers of At, keeping only terms up to the first order. From the expression 
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for p(k, 6>, n), Eq. (4.2), we have 

n+l 



p(k, 9,n+l) = exp<( i(ku n+1 + 9w n+1 ) - - 



k 2 9(n+l,i)g(n+l,j)(j>(i,j) At 2 

n+l 

+k9 [g(n+l,i)h(n+l,j) + h(n+l,i)g(n+l,j)}(j)(i,j) At 2 

i,3=0 

n+l 1 ^ 

+9 2 h(n+l,i)h(n+l,j)cf>(i,j) At 2 I. (B.l) 
i,j=o -I J 

Let us define dig as the partial derivative of g with respect to its first variable, i. e., 



so that dig(n,i) means 



dg(t,t'), 

<9i#(n, i) = — — — t =nAt ■ {B.3) 
ot t'=iAt 



A similar convention holds for the function h(t,t'). With this convention, and for the 
purpose of calculating the r.h.s of Eq. (5.1), we can approximate: 

g(n+l,i)g(n+l,j) « [g(n,i) + d ig (n,i)At] [g(n, j) + d ig (n, j)At] . (BA) 

Keeping terms only up to first order in At, we have that 

g(n + 1, i)g(n + 1, j) « g(n, i)g(n,j) + g(n, j)9i^(n, i)At + g(n, i)dig(nJ)At . (B.5a) 

g(n + l,i)h(n+lj) + h(n + l,i)g(n+ 1, j) « g(n,i)h(nj) + h(n,i)g(n,j) 

+ h{n, j)dig(n, i)At + h(n, i)d\g{n, j)At + g(n, j)dih(n, i)At + g{n, i)dih{n, j)At . 

(B.5b) 

h(n + 1, i)h(n + 1, j) h(n, i)h(n,j) + h(n, j)dih(n, i)At + /i(n, i)dih(n,j)At . (5.5c) 

du 

u n +i ^u n + -Jj-At . (B.5d) 

w n+1 ^w n + -j^At . (B.5e) 
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Hence, e. g., 

n+l n 

g(n + l,i)h(n+lJ)(P(iJ)At 2 = g(n,i)h(n,j)(P(iJ)At 2 



i,j=0 



i,j=0 
n 



+ [h(nJ)d 1 g(n,i) + g(n,i)dih(n,j)](j>(i,j)At 3 

i,3=0 
n 

+ Y 9 ( n + 1 i n )H n J) ( f>( n J) At 2 



3=0 



+ g(n, i)h(n + 1, n)4>(i, n) At z 



(5.6) 



i=0 



Note that the last term in the l.h.s sum, g(n + 1, n + l)h(n + l,n + l)<p(n + 1, n + 1) At 2 , 
is of second order in At and was therefore neglected. 
Using Eqs. (B.5), we have that to first order in At 



Ap(k, 6>, n) = p(k, 9, n) 

x j-1 + exp 

where 



at at / 2 



A= ^ 2^(n,%(n,j)^,i)At 3 + 2ai(7(n,n)^(7(?i,j)(/)(n,i) At 2 



(5.7a) 



(B.76) 



M=0 
n 



i=o 



1 

B =- Y [9i9(n,i)h(n,j) + dxhin, j)g{n,i) + dxg{n, j)h{n,i) + dih(n,i)g(n,j)]<f>(i,j) At c 



i,j=0 



+ -d-Ji{n,n) 



Y Kn, i)(f)(i, n) At 2 + ^ h(n, j)(f)(n, j) Ar 



i=0 



3=0 



^2g(n,i)(f>(i, n) At 2 + 9(n, j)(f>(n, j) At 2 



3=0 



i=0 



(B.7c) 



C=- ZdthfaiMnJMiJ) At 3 + 2d 1 h(n,n)J2Hn,j)<l>(n,j) At 2 . {Bid) 

i,j=0 j=0 

In the limit At — > 0, we can expand the exponents and keep only the leading terms. 
Thus: 



Ap(k, 6, n) = p(k, 9, n) 



. , du(n) n dw(n) 
i ( k — + 9 



dt 



dt 



- (Ak 2 + 2Bk9 + C9 2 ) 



(5.8) 
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where A, B, C are defined in Eqs. (B.7b)-(B.7d). 

In the limit n — > oo, At — > 0, all sums become integrals. Taking for example the 
coefficient of k 2 , we have: 



^ 2<9i<7(n, 1)0(71, j)(f>(i, j) At 3 + 2d l9 (n, n) 9(n, j)(f>{n, j) At 2 



i,j=0 



3=0 



y^dt'^dt" [2d 1 g(t,t')g{t,t") + 2d 1 g(t,t")g{t,t')]+ f dt' d 1 g(t 1 t)g(t, if)<f>{t, *')} At. 

(5.9a) 



The expression on the r.h.s turns out to be 

d_ r '* '* 



/" dt' / dt"g(t,t')g(t,t")cf>(t',t") 
o ./o 

which is just dot{t)/dt, where a(t) is 

a(t) = fdt> fdt" g(t, t')g(t, t") <f>(t', t") , 
Jo Jo 



At, 



(B.9b) 



(B.9c) 



This is identical to the definition of a(t) in Eq. (4.3b). 

Similar reasoning applied to the other terms in Eq. (B.8) leads us to finally rewrite 
it, in the limit n — > oo, At — > 0, as 



dp(k,9,t) j .du .dw 1 

— \ % — ■ — K ~r 1 — i — (/ — 

at 1 dt dt 2 



where a, j3 and 7 are 



da 
~dt 



dt 



dt 



p(k,9,t), (5.10) 



a(t) = / dt' / V o(t, t')<7(t, O t") , 
Jo Jo 

/3(t) = \ fdt' fdt" [g(t,t')h(t,t") + h(t,t')g(t,t")} WX), 
1 Jo Jo 

7 (t)= f dt' [ dt"h{t,t')h(t,t")(l>{t',t"), 
Jo Jo 

Appendix C 

The function p{k, 6>, n) is given by 

p(k, 9, n) = e^^W{<p u . . • , <Pn)l i= M n , i)k+h(nm * t ■ 
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(5.11a) 
(5.116) 
(5.11c) 



(C.l) 



For the Markovian dichotomous noise, we have 



where 



and 



Po(yi) = ^[S(yi-a)+d( yi +a)] , 
T A t(y\y') = \ [f(At)5(y - y') + g(At)5(y + y')] , 



/(At) = 1 + exp (-2XAt) , 
g(At) = 1 - exp (-2XAt) . 

W n (y u ...,y n ) = P Q {yi) T At (y 2 \yi) ■ ■ ■ T At (y n \y n -i) 



(C.2a) 
(C.2b) 

(C.2c) 
(C.2d) 

(C.2e) 



We must calculate the Fourier transform of this expression and substitute the appropriate 
expressions for the Fourier parameters <fi. We now define the following two functions 



1 



W+(yi,...,y n ) = -5(yi - a) T At (y 2 \yi) ■ ■ -T At (y n \y n -i) . 
W~ (yi, . . . ,y n ) = -5 {yx + a)T At (y 2 \yi) ■ ■ -T Af (y n |y n _i) • 

According to the definitions of Po(yi) and W n in Eq. (C.2), we see that 



(C.3a) 
{CM) 



Wn = W+ + W~ . 



(CA) 



Going over to the Fourier transform of the various functions, we have, using Eq. (C.2b), 



dy 



i • • • 



dy n exp [i{yi<pi H Vy n Vn)\ 



x t;%i - a) 



^fS(y 2 -yi) + ^gS(y 2 + yi] 



x T At (y 3 \y 2 ) ■ • -T At (?/ n |?/ n _i) . 



(C.5) 



And a similar equation for W n . After performing the integration on yi, we have: 



-1 poo poo 

W+(<p u ...,<p n ) = -e ia ^ J_ dy 



2 ' • • 

oo J —oo 



dy n exp [i{y 2 ip 2 H h y n <Pn)} 



^fS(y 2 -a) + ^gS(y 2 + a) 
x T At (y 3 |y 2 )---T At (y n |y n _i). (C.6a) 
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1 , • • • , <Pn 



-ia(pi 



/OO roo 
dy2 ■ ■ ■ I dy n exp [i(y 2 ¥>2 H h yn¥>n)] 
-00 J — oo 



^fS(y 2 + a) + ^gS(y 2 -a) 
x T At (y 3 |y 2 ) • • -Ta^^I^!) . (C.66) 



Using the definitions of VF + and , Eq. (C.3), we have: 



W+fa, ...,<p n ) = l e ia ^ {/Wii(V2, ••.,¥>») + </ ^"-1(^2, ...,¥>„)}• (C.7a) 
. . • , <p n ) = \ e-*** 1 {fW^fa, • • • » V») + <7 ^+-1^2, ...,¥»„)}. (C.76) 



with 



1 , 



wr = - e~ iaipi 
1 2 



We now have a 



(C.7c) 
(C.7d) 



Lemma : W+fa, ...,<p n ) = W n (-pi, . . . , -p n ). 
The proof follows immediately from (C.7) by induction on n. 
The first terms in the series j j 



are: 



w + = 1 ia^i 
1 2 



j e / ia(<fii+tfi 2 ) _|_ ^ e io(<pi-v? 2 ) 



w + = ( - 



2 ia(tp 1 +ip2+<Ps) 



+ fge 



ia(<p 1 +(p 2 -<P3) 



(C.8a) 
(C.86) 



(C.8c) 



This suggests that we can write: 



— n / 1 \ n 

k=l ^ ' 



(<Pl> •••»<?»») 



(C.9) 
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with F£ yet to be determined. Note that W~((pi, . . . , (p n ) follows from W+ . . . , (p n ) 
by virtue of the above lemma. 

Substituting the form (C.9) into the Eqs. (C.7) we find: 



k=i 



til) r-'-V^nfa, ■ ■ ■ ,vn) =\fe^ E Q) r _i - V- 1 *?- 1 ^, ...,<p n ) 

+ ??e iaVl E ( 9 J /"^'"V-^-'C-^, • • • , "¥>„) 

(C.10) 



i=i 



The first term on the r.h.s. can be rewritten as 

n-l 



n— i / -i \ n 



(C.lla) 



where the index j has been renamed fc. The second term on the r.h.s of Eq. (C.10) can be 
rewritten as 

n , x n 

^£(2) r^v-^ri 1 !-^,...,-^), (cub) 

where the new summation index /c is defined as k = j + 1. 

Comparing the two sides of Eq. (C.10) term by term with the help of Eqs. (C.ll), we 
see that we must have: 

'e^F^fe,...,^) k=l 
F£( Vl , . . . , <p n ) = I e ia ^ [F»-\<p 2 , . . . , <p n ) + F£zl(-V2, • • • , -<f n )] Kk<n 
.e^Ci 1 ^.-.-^ k = n. 

The three cases in Eq. (C.12) can be summed up as 



(C.12) 



^(^,...^ n )=e^ [Fr 1 (^ ) ... ! ^) + ^; 1 1 H 2) ...,-^)] , (C.13a) 
if we add the conventions 

F o m (^,...,^ n ) = 0. (C.136) 

F- +1 (^,...,^ n ) = 0. (C.13c) 

Finally, we note that from = ie 1 "^ 1 , we have that 



F 1 1 ( v? i)=e^ 1 



(C.14) 
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Eqs. (C.13) and (C.14) fully determine the coefficients Fg((pi, . . . , <p n ). Eqs. (C.l), 
(C.4), (C.9) and the lemma now allow us to write the function p{k, 6>, n) as 

p(k, 9, n) = exp [i(ku n + 6w n )} J2 Q) r-*" 1 ^- 1 

x [F£{<p u ...,?„) + F£(-<p u -<p n )] \ Vi=[kgM+eh(nii)]M , (C.15) 
Appendix D 

We wish to calculate the discrete version of the first derivative, Ap(k : 6>, n), defined as 
p(k,9,n + 1) - p(k,9,n), with p(k, 9, n) defined by Eq. (C.15). Rewriting (C.15) with 
n replaced by n + 1 yields 

n+l f / -1 \ n+l 

p(fc,tf,n+l)=^H-J /""V" 1 

X [F£+ V, • • • , + • • • , ~<Pn+l)] } 

n+l ✓ 1 \ n+l 

x | e *^i [F™(^ 2 , . . . , v^n+i) + F£(-<p 2 , -tp n+1 )] 

+ e~ ia ^ . . . , -<p n+1 ) + F£(p 2 , . . . , <p n+1 )] | , (D.l) 

where we have used the recurrence relations for F£, Eqs. (C.13). In Eq. (D.l), (p^ = 
[ k g ( n + 1 — % ) + 9 h ( n + 1 — % ) ] At, in accordance with our assumption that 

g(n,i) = g(n - i). 

Remembering that F™ +1 ((p 2 , • • • , V^n+i) = (see Eq. (C.13c)), we can rewrite the first 
and third terms in the r.h.s of Eq. (D.l) as: 

n+l / 1 x n+l 

I 

k= 

n / i \ n 



n-ri / -I \ n+l 

E (2 J /""V"* t^** • • • ' Vn+l) + e-^FM-W, -<Pn+l) 

= E (2 J V ia ^ F l (V2, <p n+ i) + e- ia ^FZ(- V2 , . . . , -<p n+1 )] • 

(D.2) 

In the remaining two terms in the r.h.s of Eq. (D.l), we change the summation index from 
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k to j = k — 1. Remembering that Fq , (^ 2) • • • , ¥>n+i) = ( see Eq. (C.13b)), we have: 

n+l / -i \ n+1 

E (2 J r-y -1 ^^(-^2, . . . , - Vn+1 ) + e-^F fe -(^, . . . , ^ + i)] 

= E U r "'" V "i l e ~ iaVlF i^ • • • » Vn+i) + e^Ffi-w, . . . , -*>„+!)] . 

(D.3) 

Combining (D.2) and (D.3), we obtain: 

p(M,n+l) = E U fn-k-i g k-i 
k=i ^ ' 

+ f [e-^F^(^ 2 , . . . , ^ n+1 ) + e ia ^F£(-<p 2 , . . . , - Vn+1 )] } , 



(D.4) 



with = /cg(n — z + 1) + 6>/i(n + 1 — i). From this we subtract 

p(fc, 0, n) = E o /^"V -1 [Wi, • • • M + Vi, • • • , V„)] , (^-5) 
fc=i ^ ' 

where ip^ = kg(n — i) + 9h(n — i). Note now that (fi+i = <p\ (it is for this relation that we 
assumed that g(n, i) = g(n — i) and h(n, i) = h(n — i)). Hence, Ap can be written as: 



k=i 



cn—k—l„k—l) Tpni 



f.i 



_|_ 9_ e ~ia.Vl _ I 



f_ —iaipi _|_ 9_^iaipx _ ^ 

2 2 



(D.6) 



Since ultimately we want the limit At — > of Ap/At, we now expand the expressions in 
square brackets up to first order in At. Referring to Eq. (C.2), we have 



/(At) = 2 - 2AAt + 2A 2 At 2 + . . . 
g(At) = 2XAt - 2X 2 At 2 + ... 
exp(za</?i) = 1 + iacpi a 2( ^i + • • • 



(D.7a) 
(D.7b) 
(D.7c) 
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where the last line follows from the fact that <p\ is first order in At, because ip\ = 
[kg(n) + 9h(n)]At. 

Hence, to first order in At: 

n /-, \ n 

Ap(k, 0, n) = ia Vl £ f - f n ~ k ~ V -1 [^2, • • • , ¥>n+i) - K("V?2, • • • , -v>»+i)] • 
fc=i ^ ' 

(D.8) 

This expression is not related in any obvious way to p(k, 9, n) and therefore we cannot 
write a first order differential equation for p{k,6,n). Let us therefore look at the second 
derivative of p{k,9,n), the discrete form of which is A 2 p/At 2 where A 2 p is given by 
p(k, 9, n + 2) + p(k, 9, n) — 2p{k, 0,n + 1). The expression p(k, 9, n + 2) is calculated in 
much the same way as p{k, 9,71+1) in Eqs. (D.1)-(D.4). This time, we'll need to use the 
recurrence relation Eq. (C.13) twice. The algebra is straightforward, and we obtain finally: 

p(k,9,n+2) = J2[^) r- k - X 9 k - X 

f c ia(<Po+<Pi) _j_ / 9 c ia(y 1 -ip ) _j_ f_9_ e -ia(ip 1 +ip ) _|_ 9 e ia(ip — ipi) 



+ F£(-<p 2 ,...,-<p n ) 



c2 
I 
c2 



4 

,2 



f_ e -ia(ipo+fi) _|_ [9_ e ia ( l Po->Pi) _|_ fjlgiaifo+tfii) _|_ 9_^ia{<p\-ipo) 
4 4 4 4 



(D.9) 



where, for the sake of consistency, <pi, . . . , <p> n +i are the same as in Eq. (D.l), %. e., <pi = 
[kg(n — i + 1) + 6>/i(n — i + 1)] At, and we have defined (po = [kg(n + 1) + 9h(n + 1)] At. 

Substituting Eqs. (D.4), (D.5) and (D.9) into the expression for A 2 p and expanding 
to second order in At (using Eq. (D.7)), we finally obtain: 
A 2 p(k, 9, n) =ia [((po — ip\) — 2\ip At] 



x 



.fc=i 



cn—k — 1 „k — l 



9 [Fk(^2, Vn+l) - F£(-(p 2 , -<Pn+l] 



n , x n 

~ «Vl £ (2 J f n ' k ~ 1 3 k ~ 1 l F k(^, <Pn+l) + F£(-<P2, -<Pn+l)] ■ 

(D.10) 

Note that (po — ip\ is indeed of second order in At, as it should be, because (po — ip\ = 
{k [g(n + 1) — g(n)] + 9 [h(n + 1) — h(n)] } At, and g (n + 1) — g{n) is already of first order 
in At. Comparing Eq. (D.10) with Eqs. (D.5) and (D.8), we have that 



A 2 p(k,9,n) a 2 ip 2 



At 2 



At 2 



±p(k,9,n) + 



<Po - Pi - 2Ay? At 
(p-iAt 



Ap(k,6,n) 
At 



(D.11) 
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Denoting <p = kg(t) + 9h(t), we see that 



At 2 



[kg(n) + Oh(n)} 2 ^ cp 2 



{D. 12a) 



<Po - <pi = J_ 
cpiAt (fx 



j j(n + l) -g{n) h(n + 1) - h(n) 



At 

(fxAt 



At 



1. 



(D.12c) 



where we have used the fact that </?o 5 Vi ~^ <P when n — > oo, At — > 0. Hence, in the limit 
n — > oo, At — > 0, Eq. (D.10) finally becomes: 



dt 2 V V ^ / 9t 

where 

(p = kg(t) + 6h(t) . (D.136) 
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